#include <vtkActor.h>
#include <vtkContourFilter.h>
#include <vtkDoubleArray.h>
#include <vtkNamedColors.h>
#include <vtkPointData.h>
#include <vtkPolyDataMapper.h>
#include <vtkProperty.h>
#include <vtkRenderWindow.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkRenderer.h>
#include <vtkSmartPointer.h>
#include <vtkStructuredPoints.h>

int main(int, char*[])
{
    auto vol =
            vtkSmartPointer<vtkStructuredPoints>::New();
    int dim = 11; // 26
    vol->SetDimensions(dim, dim, dim);
    vol->SetOrigin(-0.5, -0.5, -0.5);
    auto sp = 1.0 / (dim-1);
    vol->SetSpacing(sp, sp, sp);

    auto scalars =
            vtkSmartPointer<vtkDoubleArray>::New();
    scalars->SetNumberOfComponents(1);
    scalars->SetNumberOfTuples(dim * dim * dim);
    for (auto k = 0; k < dim; k++)
    {
        auto z = -0.5 + k * sp;
        auto kOffset = k * dim * dim;
        for (auto j = 0; j < dim; j++)
        {
            auto y = -0.5 + j * sp;
            auto jOffset = j * dim;
            for (auto i = 0; i < dim; i++)
            {
                auto x = -0.5 + i * sp;
                auto s = x * x + y * y + z * z - (0.4 * 0.4);
                auto offset = i + jOffset + kOffset;
                scalars->InsertTuple(offset, &s);
            }
        }
    }
    vol->GetPointData()->SetScalars(scalars);

    int z = dim/2;
    cout << "data at slice: " << z+1 << endl;
    for (int y = 0; y < dim; y++)
    {
        for (int x = 0; x < dim; x++)
        {
            double* pixel = static_cast<double*>(vol->GetScalarPointer(x,y,z));
            cout << setw(8) << pixel[0];
        }
        cout << endl;
    }
    cout << "============" << endl;

    auto contour =
            vtkSmartPointer<vtkContourFilter>::New();
    contour->SetInputData(vol);
    contour->SetValue(0, -0.05);
    contour->SetValue(1, 0.4);

    auto colors =
            vtkSmartPointer<vtkNamedColors>::New();

    auto renderer =
            vtkSmartPointer<vtkRenderer>::New();
    auto renWin =
            vtkSmartPointer<vtkRenderWindow>::New();
    renWin->AddRenderer(renderer);
    auto iren =
            vtkSmartPointer<vtkRenderWindowInteractor>::New();
    iren->SetRenderWindow(renWin);
    auto volMapper =
            vtkSmartPointer<vtkPolyDataMapper>::New();
    volMapper->SetInputConnection(contour->GetOutputPort());
    volMapper->ScalarVisibilityOff();
    auto volActor =
            vtkSmartPointer<vtkActor>::New();
    volActor->SetMapper(volMapper);
    volActor->GetProperty()->EdgeVisibilityOn();
    volActor->GetProperty()->SetColor(colors->GetColor3d("Salmon").GetData());
    renderer->AddActor(volActor);
    renderer->SetBackground(colors->GetColor3d("SlateGray").GetData());
    renWin->SetSize(512, 512);

    // interact with data
    renWin->Render();

    iren->Start();

    return EXIT_SUCCESS;
}
